###############################################################
##### Replication files for 

#######

##Run on R version 3.6.1
##NPC veriosn 1.1.0

dta<-read.csv("waves_for_npc_replication.csv")
head(dta)

dta$F3r<-1-dta$F3
dta$F4r<-1-dta$F4
dta$F5r<-1-dta$F5
dta$F7r<-1-dta$F7
dta$F14r<-1-dta$F14
dta$legal.rev<-1-dta$legal

dta$F6r<-1-dta$F6
dta$F8r<-1-dta$F8
dta$F15r<-1-dta$F15
dta$borderenr<-1-dta$borderen

lowknow<-subset(dta, dta$know_index==0)

library(NPC)

#########################################################
###### DIFFERENCE IN MEANS ###############################
#########################################################

####Low knowledge sample
##Simple herd on information
npc.simpleherd.info <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                           y.vars=c("F11", "F12" , "F13"), n.perms=1000,
                           alternative = "greater", seed=123456, 
                           comb.fun="NormalCF",test.statistic="DiffMeans", 
                           return.matrix = TRUE, FWE.adj=TRUE)
round(npc.simpleherd.info$p.value, 3)

### Bayesian herd on information
npc.bayes.info <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                      y.vars=c("F11", "F12" , "F13", "F3",
                               "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                               "violence", "goods", "turkey", "trust"), n.perms=1000,
                      alternative = "greater", seed=123456, 
                      comb.fun="NormalCF",test.statistic="DiffMeans", 
                      return.matrix = TRUE, FWE.adj=FALSE)
round(npc.bayes.info$p.value, 3)
npc.bayes.info.adjpvals<-FWE(npc.bayes.info$p.matrix, stepdown = FALSE, cfun="MinimumCF")

##Simple herd sympathetic
npc.simpleherd.symp <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                           y.vars=c("F11", "F12" , "F13"), n.perms=1000,
                           alternative = "greater", seed=123456, 
                           comb.fun="NormalCF",test.statistic="DiffMeans", 
                           return.matrix = TRUE, FWE.adj=TRUE)
round(npc.simpleherd.symp$p.value, 3)


#Bayesian herd sympathetic
npc.bayes.symp <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                      y.vars=c("F11", "F12" , "F13", "F3",
                               "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                               "violence", "goods", "turkey", "trust"), n.perms=1000,
                      alternative = "greater", seed=123456, 
                      comb.fun="NormalCF",test.statistic="DiffMeans", 
                      return.matrix = TRUE, FWE.adj=FALSE)
round(npc.bayes.symp$p.value, 3)
npc.bayes.symp.adjpvals<-FWE(npc.bayes.symp$p.matrix, stepdown = FALSE, cfun="MinimumCF")



###Opening
npc.opening <- NPC(data=lowknow, tr.var="opening_control", tr.label=1, 
                   y.vars=c("F11", "F12" , "F13", "F3",
                            "F4", "F5",  "F6", "F14", "F7", "F8", "F15"), n.perms=1000,
                   alternative = "greater", seed=123456, 
                   comb.fun="NormalCF",test.statistic="DiffMeans", 
                   return.matrix = TRUE, FWE.adj=TRUE)
round(npc.opening$p.value, 3)


#Hostile
npc.hostile <- NPC(data=lowknow, tr.var="hostile_control", tr.label=1, 
                   y.vars=c("F11", "F12" , "F13", "F3r",
                            "F4r", "F5r", "F6r", "F14r", "F7r", "F8r", "F15r"), n.perms=1000,
                   alternative = "greater", seed=123456, 
                   comb.fun="NormalCF",test.statistic="DiffMeans", 
                   return.matrix = TRUE, FWE.adj=TRUE)
round(npc.hostile$p.value, 5)




t.simple.info<-npc.simpleherd.info$T.matrix[1,]
p.simple.info<-npc.simpleherd.info$p.values
t.bayes.info<-npc.bayes.info$T.matrix[1,]
p.bayes.info<-npc.bayes.info.adjpvals[,1]
p.bayes.info2<-npc.bayes.info$p.values
t.simple.symp<-npc.simpleherd.symp$T.matrix[1,]
p.simple.symp<-npc.simpleherd.symp$p.values
t.bayes.symp<-npc.bayes.symp$T.matrix[1,]
p.bayes.symp<-npc.bayes.symp.adjpvals[,1]
p.bayes.symp2<-npc.bayes.symp$p.values
t.opening<-npc.opening$T.matrix[1,]
p.opening<-npc.opening$p.values
t.hostile<-npc.hostile$T.matrix[1,]
p.hostile<-npc.hostile$p.values


results <- matrix(NA,23,14)
results[,1] <- c(t.simple.info,rep(NA,20))
results[,2]<-c(p.simple.info, rep(NA,16))
results[,3]<-c(t.bayes.info, rep(NA,8))
results[,4]<-c(p.bayes.info, rep(NA,8))
results[,5]<-c(p.bayes.info2, rep(NA,7))
results[,6]<-c(t.simple.symp, rep(NA, 20))
results[,7]<-c(p.simple.symp, rep(NA,16))
results[,8]<-c(t.bayes.symp, rep(NA,8))
results[,9]<-c(p.bayes.symp, rep(NA,8))
results[,10]<-c(p.bayes.symp2, rep(NA,7))
results[,11]<-c(t.opening,rep(NA,12))
results[,12]<-c(p.opening)
results[,13]<-c(t.hostile,rep(NA,12))
results[,14]<-c(p.hostile)

write.csv(results, "results_npc_replication.csv")




##################################################################

####Repeat on full sample
##Simple herd on information
npc.simpleherd.info.f <- NPC(data=dta, tr.var="info_control", tr.label=1, 
                             y.vars=c("F11", "F12" , "F13"), n.perms=1000,
                             alternative = "greater", seed=123456, 
                             comb.fun="NormalCF",test.statistic="DiffMeans", 
                             return.matrix = TRUE, FWE.adj=TRUE)
round(npc.simpleherd.info.f$p.value, 3)

npc.bayes.info.f <- NPC(data=dta, tr.var="info_control", tr.label=1, 
                        y.vars=c("F11", "F12" , "F13", "F3",
                                 "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                                 "violence", "goods", "turkey", "trust"), n.perms=1000,
                        alternative = "greater", seed=123456, 
                        comb.fun="NormalCF",test.statistic="DiffMeans", 
                        return.matrix = TRUE, FWE.adj=FALSE)
round(npc.bayes.info.f$p.value, 3)
npc.bayes.info.f.adjpvals<-FWE(npc.bayes.info.f$p.matrix, stepdown = FALSE, cfun="MinimumCF")


npc.simpleherd.symp.f <- NPC(data=dta, tr.var="symp_control", tr.label=1, 
                             y.vars=c("F11", "F12" , "F13"), n.perms=1000,
                             alternative = "greater", seed=123456, 
                             comb.fun="NormalCF",test.statistic="DiffMeans", 
                             return.matrix = TRUE, FWE.adj=TRUE)
round(npc.simpleherd.symp.f$p.value, 2)

npc.bayes.symp.f <- NPC(data=dta, tr.var="symp_control", tr.label=1, 
                        y.vars=c("F11", "F12" , "F13", "F3",
                                 "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                                 "violence", "goods", "turkey", "trust"), n.perms=1000,
                        alternative = "greater", seed=123456, 
                        comb.fun="NormalCF",test.statistic="DiffMeans", 
                        return.matrix = TRUE, FWE.adj=FALSE)
round(npc.bayes.symp.f$p.value, 2)
npc.bayes.symp.f.adjpvals<-FWE(npc.bayes.symp.f$p.matrix, stepdown = FALSE, cfun="MinimumCF")


###Opening
npc.opening.f <- NPC(data=dta, tr.var="opening_control", tr.label=1, 
                     y.vars=c("F11", "F12" , "F13", "F3",
                              "F4", "F5",  "F6", "F14", "F7", "F8", "F15"), n.perms=1000,
                     alternative = "greater", seed=123456, 
                     comb.fun="NormalCF",test.statistic="DiffMeans", 
                     return.matrix = TRUE, FWE.adj=TRUE)
round(npc.opening.f$p.value, 3)


#Hostile
npc.hostile.f <- NPC(data=dta, tr.var="hostile_control", tr.label=1, 
                     y.vars=c("F11", "F12" , "F13", "F3r",
                              "F4r", "F5r", "F6r", "F14r", "F7r", "F8r", "F15r"), n.perms=1000,
                     alternative = "greater", seed=123456, 
                     comb.fun="NormalCF",test.statistic="DiffMeans", 
                     return.matrix = TRUE, FWE.adj=TRUE)
round(npc.hostile.f$p.value, 3)

npc.hostile.life.f <- NPC(data=dta, tr.var="hostile_control", tr.label=1, 
                          y.vars=c("F11", "F12" , "F13", "F3r",
                                   "F4r", "F5r", "F6r", "F14r", "F7r", 
                                   "F8r", "F15r", "knowanything_europe", "F28", 
                                   "F29", "F30"), n.perms=1000,
                          alternative = "greater", seed=123456, 
                          comb.fun="NormalCF",test.statistic="DiffMeans", 
                          return.matrix = TRUE, FWE.adj=TRUE)
round(npc.hostile.life.f$p.value, 3)


t.simple.info.f<-npc.simpleherd.info.f$T.matrix[1,]
p.simple.info.f<-npc.simpleherd.info.f$p.values
t.bayes.info.f<-npc.bayes.info.f$T.matrix[1,]
p.bayes.info.f<-npc.bayes.info.f.adjpvals[,1]
p.bayes.info.f2<-npc.bayes.info.f$p.values
t.simple.symp.f<-npc.simpleherd.symp.f$T.matrix[1,]
p.simple.symp.f<-npc.simpleherd.symp.f$p.values
t.bayes.symp.f<-npc.bayes.symp.f$T.matrix[1,]
p.bayes.symp.f<-npc.bayes.symp.f.adjpvals[,1]
p.bayes.symp.f2<-npc.bayes.symp.f$p.values
t.opening.f<-npc.opening.f$T.matrix[1,]
p.opening.f<-npc.opening.f$p.values
t.hostile.f<-npc.hostile.f$T.matrix[1,]
p.hostile.f<-npc.hostile.f$p.values


results.f <- matrix(NA,23,14)
results.f[,1] <- c(t.simple.info.f,rep(NA,20))
results.f[,2]<-c(p.simple.info.f, rep(NA,16))
results.f[,3]<-c(t.bayes.info.f, rep(NA,8))
results.f[,4]<-c(p.bayes.info.f, rep(NA,8))
results.f[,5]<-c(p.bayes.info.f2, rep(NA,7))
results.f[,6]<-c(t.simple.symp.f, rep(NA, 20))
results.f[,7]<-c(p.simple.symp.f, rep(NA,16))
results.f[,8]<-c(t.bayes.symp.f, rep(NA,8))
results.f[,9]<-c(p.bayes.symp.f, rep(NA,8))
results.f[,10]<-c(p.bayes.symp.f2, rep(NA,7))
results.f[,11]<-c(t.opening.f,rep(NA,12))
results.f[,12]<-c(p.opening.f)
results.f[,13]<-c(t.hostile.f,rep(NA,12))
results.f[,14]<-c(p.hostile.f)

write.csv(results.f, "resultsfull_npc_replication.csv")



#############################################
### On the created indices
##Simple herd on information

npc.simpleherd.info.index <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                                 y.vars=c("eusoon", "borderen"), n.perms=1000, ##Extra index included to get the pvalue
                                 alternative = "greater", seed=123456, 
                                 comb.fun="NormalCF",test.statistic="DiffMeans", 
                                 return.matrix = TRUE, FWE.adj=TRUE)
round(npc.simpleherd.info.index$p.value, 3)

### Bayesian herd on information
npc.bayes.info.index <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                            y.vars=c("eusoon", "legal", "borderen",
                                     "home"), n.perms=1000,
                            alternative = "greater", seed=123456, 
                            comb.fun="NormalCF",test.statistic="DiffMeans", 
                            return.matrix = TRUE, FWE.adj=TRUE)
round(npc.bayes.info.index$p.value, 3)

##Simple herd sympathetic
npc.simpleherd.symp.index <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                                 y.vars=c( "eusoon", "borderen"), n.perms=1000, ##Extra index included to get the pvalue
                                 alternative = "greater", seed=123456, 
                                 comb.fun="NormalCF",test.statistic="DiffMeans", 
                                 return.matrix = TRUE, FWE.adj=TRUE)
round(npc.simpleherd.symp.index$p.value, 3)


#Bayesian herd sympathetic
npc.bayes.symp.index <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                            y.vars=c("eusoon", "legal", "borderen",
                                     "home"), n.perms=1000,
                            alternative = "greater", seed=123456, 
                            comb.fun="NormalCF",test.statistic="DiffMeans", 
                            return.matrix = TRUE, FWE.adj=TRUE)
round(npc.bayes.symp.index$p.value, 3)


###Opening
npc.opening.index <- NPC(data=lowknow, tr.var="opening_control", tr.label=1, 
                         y.vars=c("eusoon", "legal", "borderen"), n.perms=1000,
                         alternative = "greater", seed=123456, 
                         comb.fun="NormalCF",test.statistic="DiffMeans", 
                         return.matrix = TRUE, FWE.adj=TRUE)
round(npc.opening.index$p.value, 3)


#Hostile
npc.hostile.index <- NPC(data=lowknow, tr.var="hostile_control", tr.label=1, 
                         y.vars=c( "eusoon", "legal.rev", "borderenr"), n.perms=1000,
                         alternative = "greater", seed=123456, 
                         comb.fun="NormalCF",test.statistic="DiffMeans", 
                         return.matrix = TRUE, FWE.adj=TRUE)
round(npc.hostile.index$p.value, 3)



t.simple.info.index<-npc.simpleherd.info.index$T.matrix[1,]
p.simple.info.index<-npc.simpleherd.info.index$p.values
t.bayes.info.index<-npc.bayes.info.index$T.matrix[1,]
p.bayes.info.index<-npc.bayes.info.index$p.values
t.simple.symp.index<-npc.simpleherd.symp.index$T.matrix[1,]
p.simple.symp.index<-npc.simpleherd.symp.index$p.values
t.bayes.symp.index<-npc.bayes.symp.index$T.matrix[1,]
p.bayes.symp.index<-npc.bayes.symp.index$p.values
t.opening.index<-npc.opening.index$T.matrix[1,]
p.opening.index<-npc.opening.index$p.values
t.hostile.index<-npc.hostile.index$T.matrix[1,]
p.hostile.index<-npc.hostile.index$p.values


results.index <- matrix(NA,9,12)
results.index[,1] <- c(t.simple.info.index,rep(NA,7))
results.index[,2]<-c(p.simple.info.index, rep(NA,4))
results.index[,3]<-c(t.bayes.info.index, rep(NA,5))
results.index[,4]<-c(p.bayes.info.index)
results.index[,5]<-c(t.simple.symp.index, rep(NA, 7))
results.index[,6]<-c(p.simple.symp.index, rep(NA,4))
results.index[,7]<-c(t.bayes.symp.index, rep(NA,5))
results.index[,8]<-c(p.bayes.symp.index)
results.index[,9]<-c(t.opening.index,rep(NA,6))
results.index[,10]<-c(p.opening.index,rep(NA,2))
results.index[,11]<-c(t.hostile.index,rep(NA,6))
results.index[,12]<-c(p.hostile.index, rep(NA, 2))

write.csv(results.index, "resultsindex_npc_replication.csv")


################################################
#### With indices from PCA

##Simple herd on information
npc.simpleherd.info.pca <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                               y.vars=c( "beeusoon", "borderenforce"), n.perms=1000, ##Extra index included to make it run
                               alternative = "greater", seed=123456, 
                               comb.fun="NormalCF",test.statistic="DiffMeans", 
                               return.matrix = TRUE, FWE.adj=TRUE)
round(npc.simpleherd.info.pca$p.value, 3)

### Bayesian herd on information
npc.bayes.info.pca <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                          y.vars=c("beeusoon", "staywork", "borderenforce", 
                                   "condhometrans", "trust2"), n.perms=1000,
                          alternative = "greater", seed=123456, 
                          comb.fun="NormalCF",test.statistic="DiffMeans", 
                          return.matrix = TRUE, FWE.adj=TRUE)
round(npc.bayes.info.pca$p.value, 3)

##Simple herd sympathetic
npc.simpleherd.symp.pca <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                               y.vars=c("beeusoon", "borderenforce"), n.perms=1000, ##Extra index included to make it run
                               alternative = "greater", seed=123456, 
                               comb.fun="NormalCF",test.statistic="DiffMeans", 
                               return.matrix = TRUE, FWE.adj=TRUE)
round(npc.simpleherd.symp.pca$p.value, 3)


#Bayesian herd sympathetic
npc.bayes.symp.pca <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                          y.vars=c("beeusoon", "staywork", "borderenforce", 
                                   "condhometrans", "trust2"), n.perms=1000,
                          alternative = "greater", seed=123456, 
                          comb.fun="NormalCF",test.statistic="DiffMeans", 
                          return.matrix = TRUE, FWE.adj=TRUE)
round(npc.bayes.symp.pca$p.value, 3)


###Opening
npc.opening.pca <- NPC(data=lowknow, tr.var="opening_control", tr.label=1, 
                       y.vars=c( "beeusoon", "staywork", "borderenforce"), n.perms=1000,
                       alternative = "greater", seed=123456, 
                       comb.fun="NormalCF",test.statistic="DiffMeans", 
                       return.matrix = TRUE, FWE.adj=TRUE)
round(npc.opening.pca$p.value, 3)


#Hostile
npc.hostile.pca <- NPC(data=lowknow, tr.var="hostile_control", tr.label=1, 
                       y.vars=c( "beeusoon",  "staywork", "borderenforce", "friends"), n.perms=1000,
                       alternative = "greater", seed=123456, 
                       comb.fun="NormalCF",test.statistic="DiffMeans", 
                       return.matrix = TRUE, FWE.adj=TRUE)
round(npc.hostile.pca$p.value, 3)



t.simple.info.pca<-npc.simpleherd.info.pca$T.matrix[1,]
p.simple.info.pca<-npc.simpleherd.info.pca$p.values
t.bayes.info.pca<-npc.bayes.info.pca$T.matrix[1,]
p.bayes.info.pca<-npc.bayes.info.pca$p.values
t.simple.symp.pca<-npc.simpleherd.symp.pca$T.matrix[1,]
p.simple.symp.pca<-npc.simpleherd.symp.pca$p.values
t.bayes.symp.pca<-npc.bayes.symp.pca$T.matrix[1,]
p.bayes.symp.pca<-npc.bayes.symp.pca$p.values
t.opening.pca<-npc.opening.pca$T.matrix[1,]
p.opening.pca<-npc.opening.pca$p.values
t.hostile.pca<-npc.hostile.pca$T.matrix[1,]
p.hostile.pca<-npc.hostile.pca$p.values


results.pca <- matrix(NA,11,12)
results.pca[,1] <- c(t.simple.info.pca,rep(NA,9))
results.pca[,2]<-c(p.simple.info.pca, rep(NA,6))
results.pca[,3]<-c(t.bayes.info.pca, rep(NA,6))
results.pca[,4]<-c(p.bayes.info.pca)
results.pca[,5]<-c(t.simple.symp.pca, rep(NA, 9))
results.pca[,6]<-c(p.simple.symp.pca, rep(NA,6))
results.pca[,7]<-c(t.bayes.symp.pca, rep(NA,6))
results.pca[,8]<-c(p.bayes.symp.pca)
results.pca[,9]<-c(t.opening.pca,rep(NA,8))
results.pca[,10]<-c(p.opening.pca,rep(NA,4))
results.pca[,11]<-c(t.hostile.pca,rep(NA,7))
results.pca[,12]<-c(p.hostile.pca, rep(NA, 2))

write.csv(results.pca, "resultspca_npc_replication.csv")

#########################################################
###### DIFFERENCE IN SUMS ###############################
#########################################################

######NPC
####Low knowledge sample
##Simple herd on information
npc.sum.simpleherd.info <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                               y.vars=c("F11", "F12" , "F13"), n.perms=1000,
                               alternative = "greater", seed=123456, 
                               comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                               return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.simpleherd.info$p.value, 3)


##Simple herd sympathetic
npc.sum.simpleherd.symp <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                               y.vars=c("F11", "F12" , "F13"), n.perms=1000,
                               alternative = "greater", seed=123456, 
                               comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                               return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.simpleherd.symp$p.value, 3)


### Bayesian herd on information
npc.sum.bayes.info <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                          y.vars=c("F11", "F12" , "F13", "F3",
                                   "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                                   "violence", "goods", "turkey", "trust"), n.perms=1000,
                          alternative = "greater", seed=123456, 
                          comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                          return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.bayes.info$p.value, 3)
npc.sum.bayes.info.adjpvals<-FWE(npc.sum.bayes.info$p.matrix, stepdown = FALSE, cfun="MinimumCF")


npc.sum.bayes.info.mcar <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                               y.vars=c("F11", "F12" , "F13", "F3",
                                        "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                                        "violence", "goods", "turkey", "trust"), n.perms=1000,
                               alternative = "greater", seed=123456, 
                               comb.fun="NormalCF",test.statistic="DiffSumObs", 
                               return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.bayes.info.mcar$p.value, 6)



#Bayesian herd sympathetic
npc.sum.bayes.symp <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                          y.vars=c("F11", "F12" , "F13", "F3",
                                   "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                                   "violence", "goods", "turkey", "trust"), n.perms=1000,
                          alternative = "greater", seed=123456, 
                          comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                          return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.bayes.symp$p.value, 3)
npc.sum.bayes.symp.adjpvals<-FWE(npc.sum.bayes.symp$p.matrix, stepdown = FALSE, cfun="MinimumCF")


npc.sum.bayes.symp.mcar <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                               y.vars=c("F11", "F12" , "F13", "F3",
                                        "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                                        "violence", "goods", "turkey", "trust"), n.perms=1000,
                               alternative = "greater", seed=123456, 
                               comb.fun="NormalCF",test.statistic="DiffSumObs", 
                               return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.bayes.symp.mcar$p.value, 8)


####Bayes both
bayes<-read.csv("waves_data_npc_bayes_combined.csv")
head(bayes)

npc.sum.bayes.symp.both <- NPC(data=bayes, tr.var="bayestreatment", tr.label=1, 
                               y.vars=c("violence_1", "goods_1", "turkey_1", "trust_1", "F6_1",
                                        "F8_1", "F15_1", "F11_1", "F12_1" , "F13_1", "F3_1",
                                        "F4_1", "F5_1", "F7_1", "F14_1","violence_2", "goods_2",
                                        "turkey_2", "trust_2", "F6_2",
                                        "F8_2", "F15_2", "F11_2", "F12_2" , "F13_2", "F3_2",
                                        "F4_2", "F5_2", "F7_2", "F14_2"), n.perms=1000,
                               alternative = "greater", seed=123456, 
                               comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                               return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.bayes.symp.both$p.value, 3)
npc.sum.bayes.symp.adjpvals.both<-FWE(npc.sum.bayes.symp.both$p.matrix, stepdown = FALSE, cfun="MinimumCF")



###Opening
npc.sum.opening <- NPC(data=lowknow, tr.var="opening_control", tr.label=1, 
                       y.vars=c( "F11", "F12" , "F13", "F3",
                                 "F4", "F5",  "F6", "F14", "F7", "F8", "F15"), n.perms=1000,
                       alternative = "greater", seed=123456, 
                       comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                       return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.opening$p.value, 3)

npc.sum.opening.mcar <- NPC(data=lowknow, tr.var="opening_control", tr.label=1, 
                            y.vars=c( "F11", "F12" , "F13", "F3",
                                      "F4", "F5",  "F6", "F14", "F7", "F8", "F15"), n.perms=1000,
                            alternative = "greater", seed=123456, 
                            comb.fun="NormalCF",test.statistic="DiffSumObs", 
                            return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.opening.mcar$p.value, 3)


#Hostile
npc.sum.hostile <- NPC(data=lowknow, tr.var="hostile_control", tr.label=1, 
                       y.vars=c( "F11", "F12" , "F13", "F3r",
                                 "F4r", "F5r", "F6r", "F14r", "F7r", "F8r", "F15r"), n.perms=1000,
                       alternative = "greater", seed=123456, 
                       comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                       return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.hostile$p.value, 5)


npc.sum.hostile.mcar <- NPC(data=lowknow, tr.var="hostile_control", tr.label=1, 
                            y.vars=c("F11", "F12" , "F13", "F3r",
                                     "F4r", "F5r", "F6r", "F14r", "F7r", "F8r", "F15r"), n.perms=1000,
                            alternative = "greater", seed=123456, 
                            comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                            return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.hostile.mcar$p.value, 5)


t.simple.info<-npc.sum.simpleherd.info$T.matrix[1,]
p.simple.info<-npc.sum.simpleherd.info$p.values
t.bayes.info<-npc.sum.bayes.info$T.matrix[1,]
p.bayes.info<-npc.sum.bayes.info.adjpvals[,1]
p.bayes.info2<-npc.sum.bayes.info$p.values
t.simple.symp<-npc.sum.simpleherd.symp$T.matrix[1,]
p.simple.symp<-npc.sum.simpleherd.symp$p.values
t.bayes.symp<-npc.sum.bayes.symp$T.matrix[1,]
p.bayes.symp<-npc.sum.bayes.symp.adjpvals[,1]
p.bayes.symp2<-npc.sum.bayes.symp$p.values
t.opening<-npc.sum.opening$T.matrix[1,]
p.opening<-npc.sum.opening$p.values
t.hostile<-npc.sum.hostile$T.matrix[1,]
p.hostile<-npc.sum.hostile$p.values


dsresults <- matrix(NA,23,14)
dsresults[,1] <- c(t.simple.info,rep(NA,20))
dsresults[,2]<-c(p.simple.info, rep(NA,16))
dsresults[,3]<-c(t.bayes.info, rep(NA,8))
dsresults[,4]<-c(p.bayes.info, rep(NA,8))
dsresults[,5]<-c(p.bayes.info2, rep(NA,7))
dsresults[,6]<-c(t.simple.symp, rep(NA, 20))
dsresults[,7]<-c(p.simple.symp, rep(NA,16))
dsresults[,8]<-c(t.bayes.symp, rep(NA,8))
dsresults[,9]<-c(p.bayes.symp, rep(NA,8))
dsresults[,10]<-c(p.bayes.symp2, rep(NA,7))
dsresults[,11]<-c(t.opening,rep(NA,12))
dsresults[,12]<-c(p.opening)
dsresults[,13]<-c(t.hostile,rep(NA,12))
dsresults[,14]<-c(p.hostile)

write.csv(dsresults, "dsresults_npc_replication.csv")


##Table A6
dsmcar<-matrix(NA, 16, 4)
dsmcar[,1]<-c(npc.sum.bayes.info.mcar$p.values)
dsmcar[,2]<-c(npc.sum.bayes.symp.mcar$p.values)
dsmcar[,3]<-c(npc.sum.opening.mcar$p.values, rep(NA,4))
dsmcar[,4]<-c(npc.sum.hostile.mcar$p.values, rep(NA,4))

write.csv(dsmcar, "ds_mcar_A6.csv")


##################################################################

####Repeat on full sample
##Simple herd on information
npc.sum.simpleherd.info.f <- NPC(data=dta, tr.var="info_control", tr.label=1, 
                                 y.vars=c("F11", "F12" , "F13"), n.perms=1000,
                                 alternative = "greater", seed=123456, 
                                 comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                                 return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.simpleherd.info.f$p.value, 3)

npc.sum.bayes.info.f <- NPC(data=dta, tr.var="info_control", tr.label=1, 
                            y.vars=c("F11", "F12" , "F13", "F3",
                                     "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                                     "violence", "goods", "turkey", "trust"), n.perms=1000,
                            alternative = "greater", seed=123456, 
                            comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                            return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.bayes.info.f$p.value, 3)
npc.sum.bayes.info.f.adjpvals<-FWE(npc.sum.bayes.info.f$p.matrix, stepdown = FALSE, cfun="MinimumCF")


npc.sum.simpleherd.symp.f <- NPC(data=dta, tr.var="symp_control", tr.label=1, 
                                 y.vars=c("F11", "F12" , "F13"), n.perms=1000,
                                 alternative = "greater", seed=123456, 
                                 comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                                 return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.simpleherd.symp.f$p.value, 2)

npc.sum.bayes.symp.f <- NPC(data=dta, tr.var="symp_control", tr.label=1, 
                            y.vars=c("F11", "F12" , "F13", "F3",
                                     "F4", "F5", "F6",  "F14", "F7","F8", "F15",
                                     "violence", "goods", "turkey", "trust"), n.perms=1000,
                            alternative = "greater", seed=123456, 
                            comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                            return.matrix = TRUE, na.rm=FALSE,  FWE.adj=FALSE)
round(npc.sum.bayes.symp.f$p.value, 2)
npc.sum.bayes.symp.f.adjpvals<-FWE(npc.sum.bayes.symp.f$p.matrix, stepdown = FALSE, cfun="MinimumCF")


###Opening
npc.sum.opening.f <- NPC(data=dta, tr.var="opening_control", tr.label=1, 
                         y.vars=c( "F11", "F12" , "F13", "F3",
                                   "F4", "F5",  "F6", "F14", "F7", "F8", "F15"), n.perms=1000,
                         alternative = "greater", seed=123456, 
                         comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                         return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.opening.f$p.value, 3)


#Hostile
npc.sum.hostile.f <- NPC(data=dta, tr.var="hostile_control", tr.label=1, 
                         y.vars=c("F11", "F12" , "F13", "F3r",
                                  "F4r", "F5r", "F6r", "F14r", "F7r", "F8r", "F15r"), n.perms=1000,
                         alternative = "greater", seed=123456, 
                         comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                         return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.hostile.f$p.value, 3)

npc.sum.hostile.life.f <- NPC(data=dta, tr.var="hostile_control", tr.label=1, 
                              y.vars=c("F11", "F12" , "F13", "F3r",
                                       "F4r", "F5r", "F6r", "F14r", "F7r", "F8r", "F15r",
                                       "knowanything_europe", "F28", 
                                       "F29", "F30"), n.perms=1000,
                              alternative = "greater", seed=123456, 
                              comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                              return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.hostile.life.f$p.value, 3)


t.simple.info.f<-npc.sum.simpleherd.info.f$T.matrix[1,]
p.simple.info.f<-npc.sum.simpleherd.info.f$p.values
t.bayes.info.f<-npc.sum.bayes.info.f$T.matrix[1,]
p.bayes.info.f<-npc.sum.bayes.info.f.adjpvals[,1]
p.bayes.info.f2<-npc.sum.bayes.info.f$p.values
t.simple.symp.f<-npc.sum.simpleherd.symp.f$T.matrix[1,]
p.simple.symp.f<-npc.sum.simpleherd.symp.f$p.values
t.bayes.symp.f<-npc.sum.bayes.symp.f$T.matrix[1,]
p.bayes.symp.f<-npc.sum.bayes.symp.f.adjpvals[,1]
p.bayes.symp.f2<-npc.sum.bayes.symp.f$p.values
t.opening.f<-npc.sum.opening.f$T.matrix[1,]
p.opening.f<-npc.sum.opening.f$p.values
t.hostile.f<-npc.sum.hostile.f$T.matrix[1,]
p.hostile.f<-npc.sum.hostile.f$p.values


dsresults.f <- matrix(NA,23,14)
dsresults.f[,1] <- c(t.simple.info.f,rep(NA,20))
dsresults.f[,2]<-c(p.simple.info.f, rep(NA,16))
dsresults.f[,3]<-c(t.bayes.info.f, rep(NA,8))
dsresults.f[,4]<-c(p.bayes.info.f, rep(NA,8))
dsresults.f[,5]<-c(p.bayes.info.f2, rep(NA,7))
dsresults.f[,6]<-c(t.simple.symp.f, rep(NA, 20))
dsresults.f[,7]<-c(p.simple.symp.f, rep(NA,16))
dsresults.f[,8]<-c(t.bayes.symp.f, rep(NA,8))
dsresults.f[,9]<-c(p.bayes.symp.f, rep(NA,8))
dsresults.f[,10]<-c(p.bayes.symp.f2, rep(NA,7))
dsresults.f[,11]<-c(t.opening.f,rep(NA,12))
dsresults.f[,12]<-c(p.opening.f)
dsresults.f[,13]<-c(t.hostile.f,rep(NA,12))
dsresults.f[,14]<-c(p.hostile.f)

write.csv(dsresults.f, "dsresultsfull_npc_replication.csv")


#############################################
### On the created indices
##Simple herd on information
npc.sum.simpleherd.info.index <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                                     y.vars=c("eusoon", "borderen"), n.perms=1000, ##Extra index included to get the pvalue
                                     alternative = "greater", seed=123456, 
                                     comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                                     return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.simpleherd.info.index$p.value, 3)

### Bayesian herd on information
npc.sum.bayes.info.index <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                                y.vars=c("eusoon","legal",  "borderen", "home"), n.perms=1000,
                                alternative = "greater", seed=123456, 
                                comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                                return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.bayes.info.index$p.value, 3)

##Simple herd sympathetic
npc.sum.simpleherd.symp.index <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                                     y.vars=c("eusoon", "borderen"), n.perms=1000, ##Extra index included to get the pvalue
                                     alternative = "greater", seed=123456, 
                                     comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                                     return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.simpleherd.symp.index$p.value, 3)


#Bayesian herd sympathetic
npc.sum.bayes.symp.index <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                                y.vars=c("eusoon","legal",  "borderen", "home"), n.perms=1000,
                                alternative = "greater", seed=123456, 
                                comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                                return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.bayes.symp.index$p.value, 3)


###Opening
npc.sum.opening.index <- NPC(data=lowknow, tr.var="opening_control", tr.label=1, 
                             y.vars=c( "eusoon", "legal", "borderen"), n.perms=1000,
                             alternative = "greater", seed=123456, 
                             comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                             return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.opening.index$p.value, 3)


#Hostile
npc.sum.hostile.index <- NPC(data=lowknow, tr.var="hostile_control", tr.label=1, 
                             y.vars=c( "eusoon", "legal.rev", "borderenr"), n.perms=1000,
                             alternative = "greater", seed=123456, 
                             comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                             return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.hostile.index$p.value, 3)



t.simple.info.index<-npc.sum.simpleherd.info.index$T.matrix[1,]
p.simple.info.index<-npc.sum.simpleherd.info.index$p.values
t.bayes.info.index<-npc.sum.bayes.info.index$T.matrix[1,]
p.bayes.info.index<-npc.sum.bayes.info.index$p.values
t.simple.symp.index<-npc.sum.simpleherd.symp.index$T.matrix[1,]
p.simple.symp.index<-npc.sum.simpleherd.symp.index$p.values
t.bayes.symp.index<-npc.sum.bayes.symp.index$T.matrix[1,]
p.bayes.symp.index<-npc.sum.bayes.symp.index$p.values
t.opening.index<-npc.sum.opening.index$T.matrix[1,]
p.opening.index<-npc.sum.opening.index$p.values
t.hostile.index<-npc.sum.hostile.index$T.matrix[1,]
p.hostile.index<-npc.sum.hostile.index$p.values


dsresults.index <- matrix(NA,9,12)
dsresults.index[,1] <- c(t.simple.info.index,rep(NA,7))
dsresults.index[,2]<-c(p.simple.info.index, rep(NA,4))
dsresults.index[,3]<-c(t.bayes.info.index, rep(NA,5))
dsresults.index[,4]<-c(p.bayes.info.index)
dsresults.index[,5]<-c(t.simple.symp.index, rep(NA, 7))
dsresults.index[,6]<-c(p.simple.symp.index, rep(NA,4))
dsresults.index[,7]<-c(t.bayes.symp.index, rep(NA,5))
dsresults.index[,8]<-c(p.bayes.symp.index)
dsresults.index[,9]<-c(t.opening.index,rep(NA,6))
dsresults.index[,10]<-c(p.opening.index,rep(NA,2))
dsresults.index[,11]<-c(t.hostile.index,rep(NA,6))
dsresults.index[,12]<-c(p.hostile.index, rep(NA, 2))

write.csv(dsresults.index, "dsresultsindex_npc_replication.csv")


################################################
#### With indices from PCA

##Simple herd on information
npc.sum.simpleherd.info.pca <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                                   y.vars=c( "beeusoon", "borderenforce"), n.perms=1000, ##Extra index included to make it run
                                   alternative = "greater", seed=123456, 
                                   comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                                   return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.simpleherd.info.pca$p.value, 3)

### Bayesian herd on information
npc.sum.bayes.info.pca <- NPC(data=lowknow, tr.var="info_control", tr.label=1, 
                              y.vars=c("beeusoon", "staywork", "borderenforce", 
                                       "condhometrans", "trust2"), n.perms=1000,
                              alternative = "greater", seed=123456, 
                              comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                              return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.bayes.info.pca$p.value, 3)

##Simple herd sympathetic
npc.sum.simpleherd.symp.pca <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                                   y.vars=c( "beeusoon", "borderenforce"), n.perms=1000, ##Extra index included to make it run
                                   alternative = "greater", seed=123456, 
                                   comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                                   return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.simpleherd.symp.pca$p.value, 3)


#Bayesian herd sympathetic
npc.sum.bayes.symp.pca <- NPC(data=lowknow, tr.var="symp_control", tr.label=1, 
                              y.vars=c("beeusoon", "staywork", "borderenforce", 
                                       "condhometrans", "trust2" ), n.perms=1000,
                              alternative = "greater", seed=123456, 
                              comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                              return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.bayes.symp.pca$p.value, 3)


###Opening
npc.sum.opening.pca <- NPC(data=lowknow, tr.var="opening_control", tr.label=1, 
                           y.vars=c("beeusoon", "staywork", "borderenforce"), n.perms=1000,
                           alternative = "greater", seed=123456, 
                           comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                           return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.opening.pca$p.value, 3)


#Hostile
npc.sum.hostile.pca <- NPC(data=lowknow, tr.var="hostile_control", tr.label=1, 
                           y.vars=c( "beeusoon",  "staywork","borderenforce",  "friends"), n.perms=1000,
                           alternative = "greater", seed=123456, 
                           comb.fun="NormalCF",test.statistic="DiffSumWithNA", 
                           return.matrix = TRUE, na.rm=FALSE,  FWE.adj=TRUE)
round(npc.sum.hostile.pca$p.value, 3)



t.simple.info.pca<-npc.sum.simpleherd.info.pca$T.matrix[1,]
p.simple.info.pca<-npc.sum.simpleherd.info.pca$p.values
t.bayes.info.pca<-npc.sum.bayes.info.pca$T.matrix[1,]
p.bayes.info.pca<-npc.sum.bayes.info.pca$p.values
t.simple.symp.pca<-npc.sum.simpleherd.symp.pca$T.matrix[1,]
p.simple.symp.pca<-npc.sum.simpleherd.symp.pca$p.values
t.bayes.symp.pca<-npc.sum.bayes.symp.pca$T.matrix[1,]
p.bayes.symp.pca<-npc.sum.bayes.symp.pca$p.values
t.opening.pca<-npc.sum.opening.pca$T.matrix[1,]
p.opening.pca<-npc.sum.opening.pca$p.values
t.hostile.pca<-npc.sum.hostile.pca$T.matrix[1,]
p.hostile.pca<-npc.sum.hostile.pca$p.values


dsresults.pca <- matrix(NA,11,12)
dsresults.pca[,1] <- c(t.simple.info.pca,rep(NA,9))
dsresults.pca[,2]<-c(p.simple.info.pca, rep(NA,6))
dsresults.pca[,3]<-c(t.bayes.info.pca, rep(NA,6))
dsresults.pca[,4]<-c(p.bayes.info.pca)
dsresults.pca[,5]<-c(t.simple.symp.pca, rep(NA, 9))
dsresults.pca[,6]<-c(p.simple.symp.pca, rep(NA,6))
dsresults.pca[,7]<-c(t.bayes.symp.pca, rep(NA,6))
dsresults.pca[,8]<-c(p.bayes.symp.pca)
dsresults.pca[,9]<-c(t.opening.pca,rep(NA,8))
dsresults.pca[,10]<-c(p.opening.pca,rep(NA,4))
dsresults.pca[,11]<-c(t.hostile.pca,rep(NA,7))
dsresults.pca[,12]<-c(p.hostile.pca, rep(NA, 2))

write.csv(dsresults.pca, "dsresultspca_npc_replication.csv")



